---
title: "Security vs. Civil Liberties: How citizens cope with threat, restriction and ideology"
author: "Felix Jäger"
date: "29.07.2022"
output:
  html_document:
    df_print: paged
---

# This file generates all results presented in the paper.

```{r set-true-if-figures-should-be-exported-as-used-in-the-document}
export_figures <- TRUE
if (export_figures & !file.exists("figures main section")){dir.create("figures main section")}
if (export_figures & !file.exists("figures appendix")){dir.create("figures appendix")}

export_numbers <- FALSE
if(export_numbers) {file.create("numbers_and_stats_in_text.txt")}
```

# Load data

```{r}
load("security_vs_civil_liberties.RData")
```

```{r packages}
library(cjoint)
library(FindIt)
library(cregg)
library(jtools)
library(sjPlot)
library(stargazer)
library(reporttools)
```

# Data preparation

```{r}
# Party vote (Sonntagsfrage)
data$party_pref[data$party_pref==""] <- NA
data$party_pref_threefold <- ifelse(data$party_pref=="AfD", "prefer_AfD",
                                    ifelse(data$party_pref=="Greens", "prefer_Greens",
                                           ifelse(is.na(data$party_pref), NA, "prefer_Other")))
data$party_pref_threefold <- factor(data$party_pref_threefold, levels = c("prefer_Other", "prefer_AfD", "prefer_Greens"))
```

```{r}
# Party affection
data$party_pref_AfD_threefold <- ifelse(data$party_pref_AfD<=3, "anti_AfD",
                                        ifelse(data$party_pref_AfD>=4 & data$party_pref_AfD<=7, "neutral_to_AfD",
                                               ifelse(data$party_pref_AfD>=8, "pro_AfD", NA)))
data$party_pref_AfD_threefold <- factor(data$party_pref_AfD_threefold, levels= c("neutral_to_AfD", "pro_AfD", "anti_AfD"))

data$party_pref_Greens_threefold <- ifelse(data$party_pref_Greens<=3, "anti_Greens",
                                        ifelse(data$party_pref_Greens>=4 & data$party_pref_Greens<=7, "neutral_to_Greens",
                                               ifelse(data$party_pref_Greens>=8, "pro_Greens", NA)))
data$party_pref_Greens_threefold <- factor(data$party_pref_Greens_threefold, levels= c("neutral_to_Greens", "pro_Greens", "anti_Greens"))
```

----------------------

# Descriptives

```{r}
col_party_three <- RColorBrewer::brewer.pal(3, "Set2")
col_party_three <- col_party_three[c(2,3,1)]
names(col_party_three) <- c("Other", "AfD", "Greens")

col_pro_anti <- RColorBrewer::brewer.pal(3, "PiYG")
col_pro_anti <- col_pro_anti[c(2,3,1)]
col_pro_anti[1] <- "#BDBDBD"
names(col_pro_anti) <- c("Neutral", "Pro", "Anti")

pch_pro_anti <- c(17, 19, 15)
```

```{r}
if(export_figures) {pdf("figures appendix/figure_desc_support_surveillance_policy.pdf", height = 4.5, width = 6)}
  barplot(table(data$vig_ter_policy)/sum(table(data$vig_ter_policy)), main="Support for surveillance policy", las=1, xlab="Support", border=F)
if(export_figures) {dev.off()}
```

```{r}
tab_survpref_party <- table(data$party_pref_threefold, data$vig_ter_policy)/rowSums(table(data$party_pref_threefold, data$vig_ter_policy))

if(export_figures) {pdf("figures appendix/figure_desc_support_surveillance_policy_by_partisanship.pdf", height = 4.5, width = 6)}
par(mar=c(4,4,2,3))
plot(0,0, pch="", xlim=c(1,10), ylim=c(0,.26), las=1, xlab = "Level of Support", ylab = "Share", bty="n", axes=F,main="Support for surveillance policy by partisanship")  
for (i in 1:3) {
  lines(1:10, tab_survpref_party[i,], col=col_party_three[i], lwd=2)
}
axis(2, las=1)
axis(1, at = c(1:10), labels = 1:10)
text(rep(10,3), tab_survpref_party[,10], labels=c("Other", "AfD", "Greens"), pos=4, col= col_party_three, xpd=T)
if(export_figures) {dev.off()}
```

```{r}
tab_party_pref <- round(sort(table(data$party_pref), decreasing = F)/sum(table(data$party_pref)),3)*100

if(export_figures) {pdf("figures appendix/figure_desc_hypothetical party vote.pdf", height = 4.5, width = 6)}
  par(mar=c(5,8,3,2))
  bp <- barplot(tab_party_pref, main="Hypothetical party vote next election", las=1, xlab="Percent", horiz = T, border=F)
  text(.8, bp, labels = tab_party_pref)
if(export_figures) {dev.off()}
```

```{r, fig.height=8, fig.width=6}
if(export_figures) {pdf("figures appendix/figure_desc_preference_afd_greens.pdf", height = 8, width = 6)}
  par(mfrow=c(2,1))
  bp <- barplot(table(data$party_pref_AfD)/sum(table(data$party_pref_AfD)), main ="Preference for AfD", las=1, border=F)
  abline(v=c((bp[3,1]+bp[4,1])/2), lty="dotted")
  abline(v=c((bp[7,1]+bp[8,1])/2), lty="dotted")
  text(c(bp[2,1],(bp[5,1]+bp[6,1])/2,bp[9,1]), c(0.4,0.4,0.4), labels = c("Anti", "Neutral", "Pro"))
  text(c(bp[2,1],(bp[5,1]+bp[6,1])/2,bp[9,1]), c(0.35,0.35,0.35), labels = (round(table(data$party_pref_AfD_threefold)/sum(table(data$party_pref_AfD_threefold)),4)*100)[c(3,1,2)])
  bp <- barplot(table(data$party_pref_Greens)/sum(table(data$party_pref_Greens)), main ="Preference for The Greens", las=1, border=F)
  abline(v=c((bp[3,1]+bp[4,1])/2), lty="dotted")
  abline(v=c((bp[7,1]+bp[8,1])/2), lty="dotted")
  text(c(bp[2,1],(bp[5,1]+bp[6,1])/2,bp[9,1]), c(0.2,0.2,0.2), labels = c("Anti", "Neutral", "Pro"))
  text(c(bp[2,1],(bp[5,1]+bp[6,1])/2,bp[9,1]), c(0.18,0.18,0.18), labels = (round(table(data$party_pref_Greens_threefold)/sum(table(data$party_pref_Greens_threefold)),4)*100)[c(3,1,2)])
if(export_figures) {dev.off()}
  
```

```{r}
round(table(data$party_pref_AfD_threefold)/sum(table(data$party_pref_AfD_threefold)),4)*100
round(table(data$party_pref_Greens_threefold)/sum(table(data$party_pref_Greens_threefold)),4)*100

if(export_numbers){write(c(paste(names(round(table(data$party_pref_Greens_threefold)/sum(table(data$party_pref_Greens_threefold)),4)*100),
                                 round(table(data$party_pref_Greens_threefold)/sum(table(data$party_pref_Greens_threefold)),4)*100),
                            paste(names(round(table(data$party_pref_AfD_threefold)/sum(table(data$party_pref_AfD_threefold)),4)*100),
                                 round(table(data$party_pref_AfD_threefold)/sum(table(data$party_pref_AfD_threefold)),4)*100)
                                ),
                         file="numbers_and_stats_in_text.txt",
                         append = T)}
```

```{r}
# Sample
table(data$age_cat)
table(data$female)
table(data$educ_cat)
data_demographics <- data.frame(data[,c("age_cat", "female", "educ_cat")])
data_demographics$age_cat <- factor(data_demographics$age_cat)
data_demographics$educ_cat <- factor(data_demographics$educ_cat, levels = c("low", "medium", "high", "Other", "Don't know"))
summary(data_demographics)
```

```{r}
colnames(data_demographics) <- c("Age", "Gender", "Education")
tableNominal(data_demographics)
```

```{r balance-table}
balance_table_list <- list(4)
for(i in 1:4) {
  data$vignette <- data[, c("vig_t_threat", "vig_t_motivation", "vig_t_target", "vig_t_party")[i] ]
  
  data$vignette <- as.factor(data$vignette)
  n_vig <- table(data$vignette)
  # Age
  a_mean <- aggregate(age ~ vignette, data = data, function(x) mean(x, na.rm=T))
  a_sd <- aggregate(age ~ vignette, data = data, function(x) sd(x, na.rm=T))
  a_n <- aggregate(age ~ vignette, data = data, length)
  a_se <- a_sd[,2]/sqrt(a_n[,2])
  
  summary(aov(age ~ vignette, data = data))
  
  # Gender
  data$gender_dummy <- data$female
  data$gender_dummy[data$gender_dummy=="divers"] <- NA
  data$gender_dummy <- as.numeric(as.factor(data$gender_dummy))
  
  g_mean <- aggregate(gender_dummy ~ vignette, data = data, function(x) mean(x, na.rm=T))
  g_sd <- aggregate(gender_dummy ~ vignette, data = data, function(x) sd(x, na.rm=T))
  g_n <- aggregate(gender_dummy ~ vignette, data = data, length)
  g_se <- g_sd[,2]/sqrt(g_n[,2])
  
  summary(aov(gender_dummy ~ vignette, data = data))
  
  # Education
  data$educ_uni_dummy <- data$educ_cat
  data$educ_uni_dummy[!data$educ_uni_dummy=="high"] <- 0
  data$educ_uni_dummy[data$educ_uni_dummy=="high"] <- 1
  data$educ_uni_dummy <- as.numeric(data$educ_uni_dummy)
  
  e_mean <- aggregate(educ_uni_dummy ~ vignette, data = data, function(x) mean(x, na.rm=T))
  e_sd <- aggregate(educ_uni_dummy ~ vignette, data = data, function(x) sd(x, na.rm=T))
  e_n <- aggregate(educ_uni_dummy ~ vignette, data = data, length)
  e_se <- e_sd[,2]/sqrt(e_n[,2])
  
  summary(aov(educ_uni_dummy ~ vignette, data = data))
  
  # Left-right
  lr_mean <- aggregate(left_right ~ vignette, data = data, function(x) mean(x, na.rm=T))
  lr_sd <- aggregate(left_right ~ vignette, data = data, function(x) sd(x, na.rm=T))
  lr_n <- aggregate(left_right ~ vignette, data = data, length)
  lr_se <- lr_sd[,2]/sqrt(lr_n[,2])
  
  summary(aov(left_right ~ vignette, data = data))
  
  
  # Combining the information into a table
  balance_table <- cbind(n_vig, a_mean, a_se, g_mean, g_se, e_mean, e_se, lr_mean, lr_se)
  balance_table <- balance_table[order(as.numeric(as.character(balance_table$vignette))),]
  balance_table <- balance_table[,c(3,2,4,5,7,8,10,11,13,14)]
  balance_table[,3:10] <- sapply(balance_table[,-c(1,2)], function (x) {round(x,3)})

  balance_table2 <- cbind(balance_table[,1:2],
                          paste0(balance_table[,3], "\n(", balance_table[,4], ")"),
                          paste0(balance_table[,5], "\n(", balance_table[,6], ")"),
                          paste0(balance_table[,7], "\n(", balance_table[,8], ")"),
                          paste0(balance_table[,9], "\n(", balance_table[,10], ")"))
  colnames(balance_table2) <- c("Vignette", "N", "Age", "Gender", "High education", "Left-Right")
  balance_table_list[[i]] <- balance_table2
}
balance_table_caption <- c("Balance check for treatment dimension threat.", 
                           "Balance check for treatment dimension terrorist motivation.",
                           "Balance check for treatment dimension scope of surveillance policy.",
                           "Balance check for treatment dimension policy proposing party.")
balance_table_label <- c("tab:balance_threat", "tab:balance_motivation", "tab:balance_policy", "tab:balance_party")
for(i in 1:4) {
  print(xtable(balance_table_list[[i]], caption=balance_table_caption[i], label=balance_table_label[i]))
}
```

```{r sample-size-interactions}
data_afd <- data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),]
data_greens <- data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),]

# H1c
print(xtable(
  table(data$vig_t_threat[!is.na(data$vig_ter_policy)], data$vig_t_target[!is.na(data$vig_ter_policy)])[c("no threat", "threat"),c("dragnet", "targeted")], caption = "Group size of values shown in Figure 2 / hypothesis H1c"))

# H2a
balance_h2a <- data.frame(cbind(
  table(data_greens$party_pref_Greens_threefold, data_greens$vig_t_party)[,"policy_from_Greens"],
  table(data_afd$party_pref_AfD_threefold, data_afd$vig_t_party)[,"policy_from_AfD"]
))
rownames(balance_h2a) <- c("neutral", "pro", "anti")
colnames(balance_h2a) <- c("Greens", "AfD")
print(xtable(balance_h2a, caption = "Group size of values shown in Figure 3 / hypothesis H2a"))

# H2b
balance_h2b <- data.frame(cbind(
  table(data_greens$party_pref_Greens_threefold, data_greens$vig_t_motivation)[,"radical_climate"],
  table(data_afd$party_pref_AfD_threefold, data_afd$vig_t_motivation)[,"right_wing"]
))
rownames(balance_h2b) <- c("Neutral", "Pro", "Anti")
colnames(balance_h2b) <- c("Greens - climate-radical", "AfD - right-wing")
print(xtable(balance_h2b, caption = "Group size of values shown in Figure 4 / hypothesis H2b"))


# H3a
balance_h3a <- data.frame(cbind(
 table(data_greens$vig_t_threat,  data_greens$vig_t_party, data_greens$party_pref_Greens_threefold)[,"policy_from_Greens","anti_Greens"],
 table(data_afd$vig_t_threat,  data_afd$vig_t_party, data_afd$party_pref_AfD_threefold)[,"policy_from_AfD","anti_AfD"]
))
rownames(balance_h3a) <- c("Control", "No pers. threat", "Pers. threat")
colnames(balance_h3a) <- c("Anti-Greens", "Anti-AfD")
print(xtable(balance_h3a, caption = "Group size of values shown in Figure 5 / hypothesis H3a"))


# H3b
balance_h3b <- data.frame(cbind(
  table(data_greens$vig_t_target,  data_greens$vig_t_party, data_greens$party_pref_Greens_threefold)[,"policy_from_Greens","pro_Greens"],
  table(data_afd$vig_t_target,  data_afd$vig_t_party, data_afd$party_pref_AfD_threefold)[,"policy_from_AfD","pro_AfD"]
))
rownames(balance_h3b) <- c("Control", "Dragnet", "Targeted")
colnames(balance_h3b) <- c("Pro Greens", "Pro AfD")
print(xtable(balance_h3b, caption = "Group size of values shown in Figure 6 / hypothesis H3b"))

```

```{r}
nrow(data)
median(as.numeric(data$Duration..in.seconds.[as.numeric(data$Duration..in.seconds.)/60<100])/60, na.rm=T)

if(export_numbers){write(c(paste("N =", nrow(data)), paste("Median response time =", median(as.numeric(data$Duration..in.seconds.[as.numeric(data$Duration..in.seconds.)/60<100])/60, na.rm=T)
)), file="numbers_and_stats_in_text.txt", append = T)}
```

----------------------

# Analysis

# Hypothesis 1

```{r}
# Hypothesis 1a+ 1b
H1ab_amce <- cjoint::amce(vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_policy),], respondent.id = "id")

# Hypothesis 1c
H1c_amce <- cjoint::amce(vig_ter_policy ~ vig_t_threat * vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_policy),], respondent.id = "id")
```

```{r}
plot(H1ab_amce)
plot(H1c_amce)
```

```{r}
H1ab_amce <- cj(data[!is.na(data$vig_ter_policy),], vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     estimate = "amce")

plot(H1ab_amce, feature_headers=T)
```

```{r}
H1ab_mm <- cj(data[!is.na(data$vig_ter_policy),], vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, estimate ="mm")
```

```{r}
if(export_figures) {pdf("figures appendix/figure_h1ab_mm.pdf")}
  plot(H1ab_mm)
if(export_figures) {dev.off()}
```

```{r}
H1ab_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_policy),])

H1c_lm <- lm(vig_ter_policy ~ vig_t_threat * vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_policy),])

plot_summs(H1ab_lm, omit.coefs = c("(Intercept)"))

if(export_figures) {pdf("figures appendix/figure_h1c_threat_target_interaction.pdf")}
  plot_model(H1c_lm , type = "pred", terms = c("vig_t_target", "vig_t_threat"))
if(export_figures) {dev.off()}

plot_model(H1ab_lm , type = "pred", terms = c("vig_t_target"))
```

```{r hypothese-1c-figure}
H1c_threat_pred <- plot_model(H1c_lm , type = "pred", terms = c("vig_t_threat", "vig_t_target"))$data
H1c_threat_pred <- data.frame(H1c_threat_pred[H1c_threat_pred$x==3,])[2:3,]

H1c_nothreat_pred <- plot_model(H1c_lm , type = "pred", terms = c("vig_t_threat", "vig_t_target"))$data
H1c_nothreat_pred <- data.frame(H1c_nothreat_pred[H1c_nothreat_pred$x==2,])[2:3,]

if(export_figures) {pdf("figures main section/figure_h1c.pdf", height = 3)}
par(mar=c(4,4,0,2))
  plot(0,0,pch="", xlim=c(1,10), ylim=c(1,2), axes=F, xlab="Predicted Values of Support", ylab="")
  for(i in 1:2) {
    points(H1c_threat_pred$predicted[i],1+i/10, col = col_pro_anti[4-i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H1c_threat_pred$conf.low[i], H1c_threat_pred$conf.high[i]), y=c(1+i/10,1+i/10), col = col_pro_anti[4-i], lwd=2)
  }
  
  for(i in 1:2) {
    points(H1c_nothreat_pred$predicted[i],2+i/10-.5, col = col_pro_anti[4-i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H1c_nothreat_pred$conf.low[i], H1c_nothreat_pred$conf.high[i]), y=c(2+i/10-.5,2+i/10-.5), col = col_pro_anti[4-i], lwd=2)
  }
  
  axis(1, at = c(1:10), labels = 1:10)
  text(x=c(0.5,0.5), y=c(1.15,2.15-.5), labels=c("Personal threat", "No personal threat"), xpd=T)
  text(x=c(rep(9.3,2)), y=c(1.1,1.2), labels=c("Dragnet", "Targeted"), xpd=T, col = col_pro_anti[3:2])
  text(x=c(rep(9.3,2)), y=c(2.1-.5,2.2-.5), labels=c("Dragnet", "Targeted"), xpd=T, col = col_pro_anti[3:2])
if(export_figures) {dev.off()}
```

```{r h1c-interaction-with-egami-framework, include = F}
h1c_fit <- CausalANOVA(vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
            int2.formula = ~ vig_t_threat : vig_t_target,
            data = data[!is.na(data$vig_ter_policy),],
            nway=2,
            family = "gaussian",
            screen=FALSE)

summary(h1c_fit)
plot(h1c_fit, type="ConditionalEffect", fac.name=c("vig_t_target","vig_t_threat"))
```

```{r, fig.height = 6.2, fig.width = 6.76}
H1ab_coef <- plot_summs(H1ab_lm, omit.coefs = c("(Intercept)"))$data

if(export_figures) {pdf("figures main section/figure_h1ab_all_coef.pdf", height = 6.2, width = 6.76)}
par(mar=c(6,9,0,0))
plot(0,0, pch ="", xlim=c(-1.5,1.3), ylim=c(1,16), axes=F, xlab="Estimate", ylab="")
  points(x=H1ab_coef$estimate[1:2], y=c(15,14), col=col_pro_anti[3], pch = 19, cex=1.5)
  for(i in 1:2) {
    lines(x=c(H1ab_coef$conf.low[i], H1ab_coef$conf.high[i]), y=c(16-i,16-i), col=col_pro_anti[3], lwd=2)
  }
  points(x=H1ab_coef$estimate[5:7], y=c(11:9), col=col_pro_anti[3], pch = 19, cex=1.5)
  for(i in 1:3) {
    lines(x=c(H1ab_coef$conf.low[i+4], H1ab_coef$conf.high[i+4]), y=c(12-i,12-i), col=col_pro_anti[3], lwd=2)
  }

  axis(2, at = c(15,14, 11, 10, 9), labels=c("Not personal", "Personal", "Islamist", "Climate-radical", "Right-wing"), las=1, tck=FALSE, lwd=0)
  axis(2, at = c(16,12), labels=c("Threat", "Terrorist motivation"), las=1, tck=FALSE, lwd=0, font=2)


  points(x=H1ab_coef$estimate[3:4], y=c(2,1), col=col_pro_anti[2], pch = 19, cex=1.5)
  for(i in 1:2) {
    lines(x=c(H1ab_coef$conf.low[i+2], H1ab_coef$conf.high[i+2]), y=c(3-i,3-i), col=col_pro_anti[2], lwd=2)
  }
  points(x=H1ab_coef$estimate[8:9], y=c(6,5), col=col_pro_anti[2], pch = 19, cex=1.5)
  for(i in 1:2) {
    lines(x=c(H1ab_coef$conf.low[i+7], H1ab_coef$conf.high[i+7]), y=c(7-i,7-i), col=col_pro_anti[2], lwd=2)
  }
  
  abline(v=0, lty="dotted")
  
  axis(1)
  axis(2, at = c(2,1, 6,5), labels=c("Dragnet", "Targeted", "AfD", "The Greens"), las=1, tck=FALSE, lwd=0)
  axis(2, at = c(3,7), labels=c("Policy scope", "Proposed by"), las=1, tck=FALSE, lwd=0, font=2)
  
  legend("bottom", legend=c("External shock", "Counter policy"), col = col_pro_anti[3:2], pch=19, xpd=T, inset = -0.25, horiz = T, bty="n")
  
if(export_figures) {dev.off()}
```

```{r}
stargazer(H1ab_lm, type="latex", no.space=T)
```

```{r}
H1ab_lm_summary <- round(coef(summary(H1ab_lm)),2)
H1ab_lm_ci <- round(confint(H1ab_lm),2)

if(export_numbers){write(c(paste("estimate policy from Greens =", H1ab_lm_summary["vig_t_partypolicy_from_Greens",1]),
                            paste("ci 2.5 policy from Greens =", H1ab_lm_ci["vig_t_partypolicy_from_Greens",1]),
                            paste("ci 97.5 policy from Greens =", H1ab_lm_ci["vig_t_partypolicy_from_Greens",2]),
                            paste("estimate policy from AfD =", H1ab_lm_summary["vig_t_partypolicy_from_AfD",1]), 
                            paste("ci 2.5 policy from AfD =", H1ab_lm_ci["vig_t_partypolicy_from_AfD",1]), 
                            paste("ci 97.5 policy from AfD =", H1ab_lm_ci["vig_t_partypolicy_from_AfD",2]), 
                            paste("estimate dragnet =", H1ab_lm_summary["vig_t_targetdragnet",1]), 
                            paste("ci 2.5 dragnet =", H1ab_lm_ci["vig_t_targetdragnet",1]), 
                            paste("ci 97.5 dragnet =", H1ab_lm_ci["vig_t_targetdragnet",2]), 
                            paste("estimate targeted =", H1ab_lm_summary["vig_t_targettargeted",1]),
                            paste("ci 2.5 targeted =", H1ab_lm_ci["vig_t_targettargeted",1]),
                            paste("ci 97.5 targeted =", H1ab_lm_ci["vig_t_targettargeted",2])
                                ),
                         file="numbers_and_stats_in_text.txt",
                         append = T)}
```


# Hypothesis 2ab

```{r}
H2b_mm <- cj(data[!is.na(data$vig_ter_policy),], vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
             by = ~ party_pref_threefold, estimate ="mm")
```

```{r}
if(export_figures) {pdf("figures appendix/figure_h2b_party_vote_mm.pdf")}
  plot(H2b_mm, group = "party_pref_threefold") + 
    aes(shape = party_pref_threefold) + # map group to `shape` aesthetic
    scale_shape_manual(values=c(15,17,19), na.translate = F) +
    scale_colour_manual(values=as.character(col_party_three[c(2,3,1)]), na.translate = F)
if(export_figures) {dev.off()}
```

```{r}
H2ab_afd_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target + party_pref_AfD_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),])
H2ab_greens_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target +  party_pref_Greens_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),])

H2b_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target +  party_pref_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_threefold),])

if(export_figures) {pdf("figures appendix/figure_h2ab_afd.pdf", height = 4.5)}
  plot_model(H2ab_afd_lm , type = "pred", terms = c("vig_t_party", "party_pref_AfD_threefold"), axis.lim = c(3,8))
if(export_figures) {dev.off()}

if(export_figures) {pdf("figures appendix/figure_h2ab_greens.pdf", height = 4.5)}
  plot_model(H2ab_greens_lm , type = "pred", terms = c("vig_t_party", "party_pref_Greens_threefold"), axis.lim = c(3,8))
if(export_figures) {dev.off()}
  
if(export_figures) {pdf("figures appendix/figure_h2b_party_vote.pdf")}
  plot_model(H2b_lm , type = "pred", terms = c("vig_t_party", "party_pref_threefold"), axis.lim = c(3,9))
if(export_figures) {dev.off()}
```

```{r}
H2ab_afd_pred <- plot_model(H2ab_afd_lm , type = "pred", terms = c("vig_t_party", "party_pref_AfD_threefold"))$data
H2ab_afd_pred <- data.frame(H2ab_afd_pred[H2ab_afd_pred$x==2,])

H2ab_greens_pred <- plot_model(H2ab_greens_lm , type = "pred", terms = c("vig_t_party", "party_pref_Greens_threefold"))$data
H2ab_greens_pred <- data.frame(H2ab_greens_pred[H2ab_greens_pred$x==3,])

if(export_figures) {pdf("figures main section/figure_h2ab.pdf", height = 3)}
par(mar=c(4,4,0,2))
  plot(0,0,pch="", xlim=c(1,10), ylim=c(1,2.5), axes=F, xlab="Predicted Values of Support", ylab="Policy proposed by")
  for(i in 1:3) {
    points(H2ab_afd_pred$predicted[i],1+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H2ab_afd_pred$conf.low[i], H2ab_afd_pred$conf.high[i]), y=c(1+i/10,1+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  for(i in 1:3) {
    points(H2ab_greens_pred$predicted[i],2+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H2ab_greens_pred$conf.low[i], H2ab_greens_pred$conf.high[i]), y=c(2+i/10,2+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  axis(1, at = c(1:10), labels = 1:10)
  text(x=c(0.5,0.5), y=c(1.2,2.2), labels=c("AfD", "Greens"), xpd=T)
  text(x=c(rep(9,3)), y=c(1.1,1.2,1.3), labels=c("Neutral", "Pro-AfD", "Anti-AfD"), xpd=T, col = col_pro_anti)
  text(x=c(rep(9,3)), y=c(2.1,2.2,2.3), labels=c("Neutral", "Pro-Greens", "Anti-Greens"), xpd=T, col = col_pro_anti)
if(export_figures) {dev.off()}
```

```{r}
H2ab_afd_mm <- cj(data[!is.na(data$vig_ter_policy),], vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
                 by = ~ party_pref_AfD_threefold, estimate ="mm")
H2ab_greens_mm <- cj(data[!is.na(data$vig_ter_policy),], vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
                    by = ~ party_pref_Greens_threefold, estimate ="mm")


if(export_figures) {pdf("figures appendix/figure_h2ab_afd_mm.pdf")}
  plot(H2ab_afd_mm , group = "party_pref_AfD_threefold")
if(export_figures) {dev.off()}

if(export_figures) {pdf("figures appendix/figure_h2ab_greens_mm.pdf")}
  plot(H2ab_greens_mm, group = "party_pref_Greens_threefold")
if(export_figures) {dev.off()}
  
```

```{r h2ab-interaction-with-egami-framework, include = F}
h2ab_fit_afd <- CausalANOVA(vig_ter_policy ~ vig_t_threat + vig_t_target + party_pref_AfD_threefold +  vig_t_party + vig_t_motivation, 
            int2.formula = ~ party_pref_AfD_threefold : vig_t_party,
            data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),],
            nway=2,
            family = "gaussian",
            screen=FALSE)

summary(h2ab_fit_afd)
plot(h2ab_fit_afd, type="ConditionalEffect", fac.name=c("party_pref_AfD_threefold","vig_t_party"))

h2ab_fit_greens <- CausalANOVA(vig_ter_policy ~ vig_t_threat + vig_t_target + party_pref_Greens_threefold +  vig_t_party + vig_t_motivation, 
            int2.formula = ~ party_pref_Greens_threefold : vig_t_party,
            data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),],
            nway=2,
            family = "gaussian",
            screen=FALSE)

summary(h2ab_fit_greens)
plot(h2ab_fit_greens, type="ConditionalEffect", fac.name=c("party_pref_Greens_threefold","vig_t_party"))
```


# Hypothesis 2c

```{r}
H2c_afd_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target + vig_t_party + party_pref_AfD_threefold * vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),])
H2c_greens_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target +  vig_t_party + party_pref_Greens_threefold * vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),])
```

```{r}
H2c_afd_pred <- plot_model(H2c_afd_lm , type = "pred", terms = c("vig_t_motivation", "party_pref_AfD_threefold"))$data
H2c_afd_pred <- data.frame(H2c_afd_pred[H2c_afd_pred$x==4,])

H2c_greens_pred <- plot_model(H2c_greens_lm , type = "pred", terms = c("vig_t_motivation", "party_pref_Greens_threefold"))$data
H2c_greens_pred <- data.frame(H2c_greens_pred[H2c_greens_pred$x==3,])

if(export_figures) {pdf("figures main section/figure_h2c.pdf", height = 3)}
par(mar=c(4,4,0,2))
  plot(0,0,pch="", xlim=c(1,10), ylim=c(1,2.5), axes=F, xlab="Predicted Values of Support", ylab="Terrorist motivation")
  for(i in 1:3) {
    points(H2c_afd_pred$predicted[i],1+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H2c_afd_pred$conf.low[i], H2c_afd_pred$conf.high[i]), y=c(1+i/10,1+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  for(i in 1:3) {
    points(H2c_greens_pred$predicted[i],2+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H2c_greens_pred$conf.low[i], H2c_greens_pred$conf.high[i]), y=c(2+i/10,2+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  axis(1, at = c(1:10), labels = 1:10)
  text(x=c(0.5,0.5), y=c(1.2,2.2), labels=c("Right-wing", "Climate-radical"), xpd=T, pos=4)
  text(x=c(rep(9,3)), y=c(1.1,1.2,1.3), labels=c("Neutral", "Pro-AfD", "Anti-AfD"), xpd=T, col = col_pro_anti)
  text(x=c(rep(9,3)), y=c(2.1,2.2,2.3), labels=c("Neutral", "Pro-Greens", "Anti-Greens"), xpd=T, col = col_pro_anti)
if(export_figures) {dev.off()}
```


# Hypothesis 3

```{r, fig.width=10}
H3a_lm <- lm(vig_ter_policy ~ vig_t_threat * party_pref_threefold *  vig_t_party + vig_t_target + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_threefold),])
H3b_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target * party_pref_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_threefold),])
```

```{r, fig.height=6}
H3a_afd_lm <- lm(vig_ter_policy ~ vig_t_threat * party_pref_AfD_threefold *  vig_t_party + vig_t_target + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),])
H3a_greens_lm <- lm(vig_ter_policy ~ vig_t_threat *  party_pref_Greens_threefold *  vig_t_party + vig_t_target + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),])
```

```{r}
H3a_afd_pred <- plot_model(H3a_afd_lm , type = "pred", terms = c("vig_t_threat", "party_pref_AfD_threefold", "vig_t_party"))$data
H3a_afd_pred <- data.frame(H3a_afd_pred[H3a_afd_pred$facet=="policy_from_AfD",])
H3a_afd_pred <- H3a_afd_pred[H3a_afd_pred$group=="anti_AfD",]

H3a_greens_pred <- plot_model(H3a_greens_lm , type = "pred", terms = c("vig_t_threat", "party_pref_Greens_threefold", "vig_t_party"))$data
H3a_greens_pred <- data.frame(H3a_greens_pred[H3a_greens_pred$facet=="policy_from_Greens",])
H3a_greens_pred <- H3a_greens_pred[H3a_greens_pred$group=="anti_Greens",]

if(export_figures) {pdf("figures main section/figure_h3a.pdf", height = 3)}
par(mar=c(4,4,0,2))
  plot(0,0,pch="", xlim=c(1,10), ylim=c(1,2.5), axes=F, xlab="Predicted Values of Support", ylab="")
  for(i in 1:3) {
    points(H3a_afd_pred$predicted[i],1+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H3a_afd_pred$conf.low[i], H3a_afd_pred$conf.high[i]), y=c(1+i/10,1+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  for(i in 1:3) {
    points(H3a_greens_pred$predicted[i],2+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H3a_greens_pred$conf.low[i], H3a_greens_pred$conf.high[i]), y=c(2+i/10,2+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  axis(1, at = c(1:10), labels = 1:10)
  text(x=c(0.5,0.5), y=c(1.2,2.2), labels=c("Anti-AfD", "Anti-Greens"), xpd=T)
  text(x=c(rep(9,3)), y=c(1.1,1.2,1.3), labels=c("Control", "No pers. threat", "Pers. threat"), xpd=T, col = col_pro_anti)
  text(x=c(rep(9,3)), y=c(2.1,2.2,2.3), labels=c("Control", "No pers. threat", "Pers. threat"), xpd=T, col = col_pro_anti)
if(export_figures) {dev.off()}
```

```{r}
round(H3a_greens_pred$predicted[3]-H3a_greens_pred$predicted[1],2)
round(H3a_afd_pred$predicted[3]-H3a_afd_pred$predicted[1],2)

if(export_numbers){write(c(paste("H3a difference Greens (pers. threat - control) = ", round(H3a_greens_pred$predicted[3]-H3a_greens_pred$predicted[1],2)),
                           paste("H3a difference AfD (pers. threat - control) = ", round(H3a_afd_pred$predicted[3]-H3a_afd_pred$predicted[1],2))
                                ),
                         file="numbers_and_stats_in_text.txt",
                         append = T)}
```


```{r, fig.height=6}
H3b_afd_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target * party_pref_AfD_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_AfD_threefold),])
H3b_greens_lm <- lm(vig_ter_policy ~ vig_t_threat + vig_t_target *  party_pref_Greens_threefold *  vig_t_party + vig_t_motivation, 
     data = data[!is.na(data$vig_ter_policy) & !is.na(data$party_pref_Greens_threefold),])
```

```{r}
H3b_afd_pred <- plot_model(H3b_afd_lm , type = "pred", terms = c("vig_t_target", "party_pref_AfD_threefold", "vig_t_party"))$data
H3b_afd_pred <- data.frame(H3b_afd_pred[H3b_afd_pred$facet=="policy_from_AfD",])
H3b_afd_pred <- H3b_afd_pred[H3b_afd_pred$group=="pro_AfD",]

H3b_greens_pred <- plot_model(H3b_greens_lm , type = "pred", terms = c("vig_t_target", "party_pref_Greens_threefold", "vig_t_party"))$data
H3b_greens_pred <- data.frame(H3b_greens_pred[H3b_greens_pred$facet=="policy_from_Greens",])
H3b_greens_pred <- H3b_greens_pred[H3b_greens_pred$group=="pro_Greens",]

if(export_figures) {pdf("figures main section/figure_h3b.pdf", height = 3)}
par(mar=c(4,4,0,2))
  plot(0,0,pch="", xlim=c(1,10), ylim=c(1,2.5), axes=F, xlab="Predicted Values of Support", ylab="")
  for(i in 1:3) {
    points(H3b_afd_pred$predicted[i],1+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H3b_afd_pred$conf.low[i], H3b_afd_pred$conf.high[i]), y=c(1+i/10,1+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  for(i in 1:3) {
    points(H3b_greens_pred$predicted[i],2+i/10, col = col_pro_anti[i], pch = pch_pro_anti[i], cex=1.5)
    lines(x=c(H3b_greens_pred$conf.low[i], H3b_greens_pred$conf.high[i]), y=c(2+i/10,2+i/10), col = col_pro_anti[i], lwd=2)
  }
  
  axis(1, at = c(1:10), labels = 1:10)
  text(x=c(0.5,0.5), y=c(1.2,2.2), labels=c("Pro-AfD", "Pro-Greens"), xpd=T)
  text(x=c(rep(10,3)), y=c(1.1,1.2,1.3), labels=c("Control", "Dragnet", "Targeted"), xpd=T, col = col_pro_anti)
  text(x=c(rep(10,3)), y=c(2.1,2.2,2.3), labels=c("Control", "Dragnet", "Targeted"), xpd=T, col = col_pro_anti)
if(export_figures) {dev.off()}
```

```{r}
round(H3b_greens_pred$predicted[2]-H3b_greens_pred$predicted[1],2)
round(H3b_afd_pred$predicted[2]-H3b_afd_pred$predicted[1],2)

if(export_numbers){write(c(paste("H3b difference Greens (dragnet - control) = ", round(H3b_greens_pred$predicted[2]-H3b_greens_pred$predicted[1],2)),
                           paste("H3b difference AfD (dragnet - control) = ", round(H3b_afd_pred$predicted[2]-H3b_afd_pred$predicted[1],2))
                                ),
                         file="numbers_and_stats_in_text.txt",
                         append = T)}
```

# Manipulation check

```{r}
policy_goal_amce <- cj(formula = vig_ter_goal ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_efficacy),], id = ~id, estimate = "amce")
concern_personal_amce <- cj(formula = vig_ter_percep1 ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_efficacy),], id = ~id, estimate = "amce")
concern_social_amce <- cj(formula = vig_ter_percep2 ~ vig_t_threat + vig_t_target + vig_t_motivation +  vig_t_party, 
     data = data[!is.na(data$vig_ter_efficacy),], id = ~id, estimate = "amce")
```

```{r}
if(export_figures) {pdf("figures appendix/manipulation_check_target.pdf")}  
  plot(policy_goal_amce) # Targeted has an effect on feeling restricted
if(export_figures) {dev.off()}
  
if(export_figures) {pdf("figures appendix/manipulation_check_threat_personal.pdf")}  
  plot(concern_personal_amce) # No effect
if(export_figures) {dev.off()}
  
if(export_figures) {pdf("figures appendix/manipulation_check_threat_societal.pdf")}  
  plot(concern_social_amce) # No effect
if(export_figures) {dev.off()}
```

